Distributions of the Baryon Fraction on Large Scales in the 

Universe 



Ping He 1,2 , Long-Long Feng 3 ' 1 , and Li-Zhi Fang 2 
ABSTRACT 

in 
o i 

The nonlinear evolution of a system consisting of collisional baryons and col- 
lisionless dark matter is generally characterized by strong shocks and disconti- 
nuities in the baryon fluid. The baryons slow down significantly at postshock 
areas of gravitational strong shocks, which can occur in high overdense as well 
as low overdense regions. On the other hand, the shocks do not affect the col- 
lapse of the dark matter. Consequently, the baryon fraction would be nonuni- 
form on large scales. We studied these phenomena with simulation samples pro- 
duced by the weighted essentially nonoscillatory (WENO) hybrid cosmological 
hydrodynamic/iV-body code, which is effective at capturing shocks and complex 
structures with high precision. We find that the baryon fraction in high mass 
[ density regions is lower on average than the cosmic baryon fraction, and many 

baryons accumulate in the regions with moderate mass density to form a high 
baryon fraction phase (HBFP). In dense regions with p/p > 100, which are the 
possible hosts for galaxy clusters, the baryon fraction can be lower than the cos- 
mic baryon fraction by about 10%-20% at z ~ 0. We also find that at z < 2, 
almost all the HBFP gas locates in the regions with mass density 0.5 < p/p < 5 
and temperature T > 10 5 K, and conversely, almost all the gas in the areas of 
0.5 < p/p < 5 and with temperature T > 10 5 K has high baryon fraction. Our 
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simulation samples show that about 3% of the cosmic baryon budget was hidden 
in the HBFP at redshift z = 3, while this percentage increases to about 14% at 
the present day. The gas in the HBFP cannot be detected either by Lya forests 
of QSO absorption spectra or by soft X-ray background. That is, the HBFP 
would be missed in the baryon budget given by current observations. 

Subject headings: cosmology: theory - intergalactic medium - large-scale struc- 
ture of the universe - methods: numerical - shock waves 
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1. Introduction 

Although the universe is dominated by dark matter and dark energy, the observed lumi- 
nous universe exists in the form of baryonic matter. The primordial nucleosynthesis predicts 
0.015 < Q b h 2 < 0.021 (Walker et al. 1991; Esposito et al. 2000) The best fitting of cosmo- 
logical parameters with the temperature fluctuations of the cosmic microwave background 
(CMB) radiation and large-scale structure clustering shows that the mass density of bary- 
onic mater is tt b = 0.0224 ± 0.0009 h~ 2 and the total matter density Vt m = 0.135l^ h~ 2 
(Bennett et al. 2003). Therefore, the cosmic baryon fraction is f c = Q b /tt m = 0.1661o;q]J. 

In the linear evolution of gravitational clustering perturbations, the density and velocity 
distributions of the baryonic gas (or intergalactic medium [IGM]) are the same as those of 
the dark matter field point-by-point on scales larger than the Jeans length. Even if the IGM 
is initially distributed differently from the dark matter, a linear growth mode will lead to 
the same distribution of the IGM as of dark matter (Bi et al. 1992; Fang et al. 1993; Nusser 
2000; Nusser & Haehnelt 1999). Thus, in the linear regime, although the density distributions 
of both IGM and dark matter are inhomogeneous, the baryon fraction is uniform on scales 
larger than the Jeans length. However, observations show that the distribution of the baryon 
fraction probably is not uniform. X-ray measurements have revealed that the baryon fraction 
in galaxy clusters is less than the prediction of primordial nucleosynthesis (Ettori & Fabian 
1999). This discrepancy is more serious in the cores of clusters (Sand et al. 2003). Even 
considering a depletion of baryons at the virial radius, the baryon fraction in galaxy clusters 
is still less than the predicted value (Ettori 2003). This indicates that in the nonlinear 
evolved fields, the baryon fraction distribution /^(x) is spatially dependent, not equal to f c 
everywhere. 

In this paper, we explore the formation and evolution for the nonuniform distribution 
of baryon fraction. On large scales, the nonuniformity of the baryon fraction is a result of 
the statistical discrepancy of baryonic gas from the underlying dark matter during nonlinear 
evolution. Recently, the statistical decoupling between baryonic gas and dark matter fields 
has been studied, both theoretically and numerically (Feng et al. 2003; He et al. 2004; 
Pando et al. 2004). These studies found that, although in both the linear and the nonlinear 
regimes the evolution of the IGM is dynamically governed by the gravity of the underlying 
dark matter field, the effects of the gravity in different regimes are very different. In the 
former, the gravity of the dark matter ensures that baryonic gas follows dark matter point- 
by-point, while in the latter, the gravity of the dark matter will inevitably lead to the 
decoupling of the baryonic gas from the dark matter on scales larger than the Jeans length. 
That is to say, the linear dynamical behavior of the baryonic matter can be simply obtained 
from the dark matter field via a similarity mapping (e.g., Kaiser 1986), while in the nonlinear 
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regime, the similarity is broken. Obviously, a direct consequence of the baryon - dark matter 
discrepancy is the inhomogeneity of the baryon fraction distribution. We follow this clue to 
study the properties of the deviation of /&(x) from f c . 

The outline of this paper is as follows. In §2 we describe the dynamical mechanism 
leading to the deviation of the baryon fraction from the cosmic value, together with the 
predictions from this mechanism. In §3 we present the samples used to numerically study 
the baryon fraction. In §4 we investigate the statistical features of the field /&(x), and 
compare them with the predictions. Finally, conclusions and discussions are given in §5. 

2. The mechanism leading to non-uniform distribution of baryon fraction 

A dynamical mechanism of separating baryonic gas from dark matter was addressed 
in the early study of structure formation (Shandarin & Zel'dovich 1989). Because the dark 
matter particles are collisionless, the velocities of the dark matter particles are multivalued 
at the intersection of the dark matter particle trajectories. On the other hand, the IGM, 
as an ideal fluid, has a single-value velocity field. Thus, discontinuities, such as shocks or 
complex structures, will develop in the density and velocity fields of gas at the intersection 
of the dark matter particle trajectories. This leads to the decoupling between the mass and 
velocity fields of the IGM and the dark matter. This feature can also be seen with the self- 
similar solution of spherical collapse under the self-gravity of baryonic gas and dark matter 
given by Bertschinger (1985). It shows that an outgoing shock is always formed during the 
infall of baryons. 

A shock is actually a common feature of IGM fluid. Although the cosmic baryonic gas 
is a Navier-Stokes fluid, the dynamic behavior of the IGM is dominated by the gravity of the 
growth modes of the dark matter. It has been recognized that the growth mode dynamics 
of cosmic baryonic gas can be approximately described by the random-force-driven Burgers 
equation (Gurbatov et al. 1989; Vergassola et al. 1994; Berera & Fang 1994; Jones 1999; 
Matarrese & Mohayaee 2002; Pando et al. 2004; He et al. 2004). The dynamical behavior 
of a Burgers fluid depends on two characteristic scales: (1) the correlation length of the 
random force (in our case, the gravity of dark matter), and (2) the Jeans length. When 
the former is larger than the latter, Burgers turbulence develops in the baryonic gas in the 
non-linear regime. That is, for the cosmic initial perturbations, which contain fluctuations 
on small as well as larger scales, the Burgers turbulence will definitely develop in nonlinear 
regime. The Burgers turbulence is qualitatively different from Navier-Stokes' turbulence. 
The latter generally consists of vortices on various scales, while the former is a collection 
of shocks (Polyakov 1995; Bouchaud et al. 1995; Yakhot 1998; Lassig 2000; Davoudi et al. 
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2001). These features arise because the velocity field is irrotational. Thus, the IGM velocity 
field in the nonlinear regime can be understood as a field consisting of shocks. 

Gas should slow down after passing through shocks, and postshock gas should have 
higher density. On the other hand, dark matter is not affected by the shocks. The nonuni- 
form distribution of baryon fraction is then caused by the postshock slowdown of gas. Con- 
sequently, Fb will be larger than 1 in the regions of postshock, which can be in high-density 
(p/p > 5) and low-density (p/p ~ 1 - 5) regions. 

Since the initial density perturbations of dark matter contain components that have 
correlation scales larger than the Jeans length at low density regions, the Burgers turbulence 
and the shocks can happen in low as well as high overdense regions. This property has also 
been shown in simulation results of He et al. (2004), who found that the shock heating is 
significant in the density regions of p/p ~ 1 — 5 at redshift z ~ 0. Therefore, the postshock 
slowdown mechanism will take place in both high and low overdense regions. 

Considering gas is moving from low- to high-density regions, the postshock slowdown 
mechanism leads to the distributions of /b(x) having the following features. 

1. Since all massive halos, at whatever redshift, formed by gravitational collapse, the 
baryon shortage in high overdense regions should be common at low as well as high redshifts. 

2. Since baryons are detained on the way from low-density to high-density (collapsed) 
regions, the /fc(x) > f c regions should be located in the low or moderate overdense areas 
surrounding high overdense area, like massive halos. 

3. Since baryons are detained in the postshock regions, the /b(x) > f c regions tend to 
be located in the higher temperature regions. 

4. Since baryons are frequently hampered by shocks, the velocities of baryonic gas are 
statistically lower than those of dark matter. 

The last point has been analyzed in Pando et al. (2004), who showed that the probability 
distribution of baryons with large peculiar velocities is much less than that of dark matter. 

3. Simulation samples 

The hydrodynamic equations of the baryonic gas in the universe is the typical Navier- 
Stokes equation (Feng et al. 2004). As discussed in the last section, the baryonic gas in 
the nonlinear regime is characterized by (1) regions with discontinuities and strong shocks 
and (2) regions with smooth and simple variations of the field between the discontinuities. 
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Therefore, an optimal simulation scheme should be effective at capturing shock and discon- 
tinuity transitions and in the meantime, at calculating piecewise smooth functions with a 
high resolution. 

For these reasons we do not use numerical schemes based on smoothed particle hydro- 
dynamic (SPH) algorithms. It is well known that one of the main challenges to the SPH 
scheme is how to handle shocks or discontinuities, because the nature of SPH is to smooth 
the fields considered (e.g., B0rve et al. 2001; Omang et al. 2003). Instead, we apply an 
Eulerian approach to simulate the IGM. Among the popular algorithms of high- resolution 
shock capturing are the total variation diminishing (TVD) scheme (Harten 1983) and the 
piecewise parabolic method (PPM; Collella & Woodward 1984). Both schemes start from 
the integral form of conservation laws of Euler equations and compute the flux vector based 
on cell averages (finite volume scheme). These methods are able to produce relatively sharp, 
nonoscillatory shock transitions. However, the TVD scheme generally degenerates to first- 
order accuracy at locations of smooth extrema (Godlewski & Raviart 1996), and this problem 
is serious in calculating the difference between hydrodynamic quantities on both sides of the 
shock when the Mach number of a gas is high. This is exactly the case that is encountered 
in the gravitationally coupled IGM and dark matter system. 

Later, the essentially nonoscillatory (ENO) and then the weighted essentially nonoscil- 
latory (WENO) schemes were proposed as an improvement over the TVD and PPM schemes 
(Harten et al. 1986; Shu 1998; Fedkiw et al. 2003; Shu 2003). It has been shown that for solv- 
ing problems governed by the high Reynolds number Navier-Stokes equations, the WENO 
is more efficient than the TVD and PPM schemes (Shi et al. 2003; Zhang et al. 2003). The 
WENO has been successfully applied to problems of astrophysical hydrodynamics, including 
stellar atmospheres (Del Zanna et al. 1998), high Reynolds number compressible flows with 
supernovae (Zhang et al. 2003), and high Mach number astrophysical jets (Carrillo et al. 
2003). In the context of cosmological applications, the WENO scheme has proved especially 
adept at handling the Burgers equation (Shu 1999). Hence, we believe that the WENO 
scheme would be effective at studying the problems sensitive to strong shocks during the 
cosmological gravitational clustering. 

Recently, a hybrid hydrodynamic/iV-body code based on the WENO scheme was devel- 
oped and has passed typical reliability tests, such as the Sedov blast wave and the formation 
of Zel'dovich pancakes (Feng et al. 2004). The code has been tested with capturing grav- 
itational shocks during the large-scale structure formation (He et al. 2004). The WENO 
algorithm on cosmological problems can be found in Feng et al. (2004). 

For the purpose of this paper, we use the same simulation samples as in He et al. (2004). 
The simulations were performed in a periodic, cubic box of size 25 /i _1 Mpc with a 192 3 grid 
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and an equal number of dark matter particles. We use the clouds-in-cell method for mass 
assignment and interpolation and adopt the seven-point finite difference to approximate the 
Laplacian operator. The simulations start at redshift z = 49, and the results are output at 
redshifts z=6, 4, 3, 2, 1, 0.5, and 0. The atomic processes, including ionization, radiative 
cooling, and heating, are modeled similarly as in Cen (1992) in a plasma of hydrogen and 
helium of primordial composition (X = 0.76, Y = 0.24). Processes such as star formation 
and feedback due to SNe and AGN activities are not taken into account yet. A uniform 
UV background of ionizing photons is assumed to have a power-law spectrum of the form 
J{v) = J 2 i x 10~ 21 (z//V_e7)~ a ergs s _1 cm~ 2 sr _1 Hz _1 , with parameters J 2 i = 1.0 and a — 1. 
The photoionizing flux is suddenly switched on at z > 10 to heat the gas and reionize the 
universe. 

For statistical studies, we randomly sampled 500 one-dimensional fields from the simu- 
lation results at each redshift. Each one-dimensional sample, of size L=25 h~ l Mpc, contains 
192 data points. For each of these points, information about the mass density and peculiar 
velocity of dark matter, as well as the mass density, peculiar velocity, and temperature of 
the IGM, is recorded. We emphasize that, since we mostly focus on one-point statistics, our 
results are dependent only on the fair sampling of these data points, which is not relevant 
to the sample dimensions. 

4. Statistical properties of the baryon fraction 

4.1. An example of the baryon fraction distribution 

Figure 1 is an example of the one-dimensional spatial distribution of Fj,(x), which is the 
normalized baryon fraction F b (x) = fb(x)/f c , and x is the comoving coordinate. We also 
show in Figure 1 the mass densities pdm and Pigm of dark matter and baryonic matter and 
the temperature T of the baryons. The mass densities p^M and Pigm are in units of their 
p DM and ~Pigmi respectively. Figure 1 clearly shows that F b (x) is significantly nonuniform. 
Some small-scale deviation of F b (x) from unity can be explained by the Jeans smoothing. 
However, there are deviations on scales of one to a few h^ 1 Mpc, which is larger than the 
corresponding Jeans length of the baryonic gas. 

In Figure 1 the mass density on the left side of the simulation box is higher than on the 
right side, and from the comparison of density distributions at z — 0.5 and 0, one can see 
that matter is infalling from the right to the left region. We see a shock at about x — 11 
h^ 1 Mpc at z — 0.5. The temperature of preshock gas is ~10 3 K and increases to ~10 7 K 
after passing through the shock. That is, the temperature increases by a factor of ~ 10 4 . 
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According to the shock theory of a polytropic gas with index 7 (Landau & Lifshitz 1959), 
the Mach number should be M ~ VlO 4 ~ 100. This value is reasonable, since under such a 
strong shock, the density of baryons is enhanced by a factor of ~ (7 + l)/(7 — 1) = 4, while 
the mass density of dark matter is not affected by the shock. Therefore, F b (x) can be as 
high as 4. The positions of this shock at time z = 0.5 and indicate that it is an outgoing 
shock with respect to the high-density area. At z — 0.5, the peak of F b (x) is just located 
on the postshock side or high-temperature side. At z — 0, more high peaks of F b (x) develop 
behind the shock. Namely, more baryonic matter is detained by the shock. This result is 
consistent with the picture of postshock slow down of the baryon flow described in §2. 

This picture can be more clearly seen in Figures 2 and 3. Figure 2 presents two- 
dimensional contours of the baryonic gas and dark matter densities. One can see at the 
lower part of the plots that a massive halo on the scale of ~5 hr 1 Mpc is formed. The 
Pigm > 1 region is obviously larger than the dark matter counterpart at puu > 1- This 
means that more baryons remain in the low overdense area. This discrepancy cannot be 
caused by the Jeans diffusion. Figure 3 gives the two-dimensional contours of temperature 
and Pigm/pdm (i-e., F b ) of the same slice in Figure 2. Figure 3 shows that the size of the 
high-temperature region is about the same as that of the Pigm > 1 region. From the right 
panel of Figure 3, we see that the contours of F b > 1 are located outside the center of the 
object. Conversely, the central part of this object, where the mass density is higher than 
that at the boundary, has only F b < 1. This feature is very common. All the clustered 
structures like the one in Figures 2 and 3 have F b < 1 in their center and are surrounded by 
F b > 1 regions. That is, baryons accumulate in low and moderate (dark matter) overdense 
areas. 



4.2. Baryon fraction - density and baryon fraction - temperature relations 

In Figure 4, we show the scatter plots between the baryon fraction and the dark matter 
density at redshifts z=6, 4, 2, 1, 0.5, and 0. Each panel in Figure 4 contains 19,200 data 
points of the randomly drawn 500 one-dimensional samples at each redshift. It is expected 
that, in each panel, most of the data points are distributed around F b ~ 1. If baryons 
underwent only linear evolution, all the points should be on the line of F b ~ 1 for all pdm- 
However, the scatters around F b = 1 are significant, among which the relatively small scatters 
can be explained by the Jeans smoothing, while the points with F b ^> 1 should be attributed 
to strong shocks or discontinuity transitions. These points are especially prominent in the 
regions of p DM < 5. 

We can see from Figure 4 that the scatter of F b at redshifts z ~ 2 is as significant as that 
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at low redshifts. Only at z > 4 does the F b scatter become small. This is consistent with 
the prediction in §2 that the deviation of F b from 1 at early times should be as prominent 
as that at later times. This feature is much different from the virialization of baryons in 
collapsed massive halos. The virialization process generally takes place during the formation 
of the halos, while the shortage of baryons in massive halos proceeds prior to the formation 
of the halos. 

Figure 5 gives the relation between the baryon fraction and the dark matter density at 
z = for baryon and dark matter fields smoothed on scales 0.26, 1.04, and 4.17 h^ 1 Mpc, 
respectively (the method of smoothing is given in Appendix A). By this treatment, the linear 
size of each cell will increase by factors of 2, 2 3 , and 2 5 , respectively. We see that the scatter 
of the smoothed F b on the scale 1.04 h" 1 Mpc remains nearly the same as that on the scale 
0.13 hr 1 Mpc even when p DM < 1. Since the size of high-F fe regions is generally larger than 
1 hr 1 Mpc (see Figs 1, 2, and 3), the scatter cannot be erased by smoothing on scales of 
about 1 hr 1 Mpc. Therefore, the scatter is intrinsic. 

As a quantitative comparison for the dependence of F b on the dark matter density, in 
Table 1 we list the mean baryon fraction in each density interval at several redshifts. We 
see from Table 1 that at all redshifts from to 4, the mean F b is always in the range of 0.7 
to 0.9 for the very high densities. Since such high-density areas are the possible sites for the 
formation of galaxy clusters, this result provides a cosmological explanation of the shortage 
of baryons in dense dark halos, like galaxy clusters. We can see that the baryon fraction of 
clusters can be as low as ~ 0.8 in general. 

In Figure 6, we show the relations between baryon fraction and temperature of baryonic 
gas at redshifts z=4, 2, 1, and 0. Each panel in Figure 6 contains 19,200 data points of the 
randomly drawn 500 samples at each redshift. The z-evolution of F b in T-space is different 
from that of p^M-space. In Figure 4 the high-i 7 }, points are always in the density range 
0.5 < pdm < 5 for all redshifts considered, while in Figure 6, the high-Fj, points are mostly 
located at 10 4 < T < 10 5 K for z = 4 and T > 10 5 K for z = 0. That is, the baryonic gas 
with F b > 1 generally tends to lie in high-temperature (T > 10 5 K) areas. This is consistent 
with the picture that strong shocks can reduce the flow of the baryons in postshock areas 
and meanwhile increase the gas temperature by a factor of 10 2 -10 4 (Fig. 1). 

4.3. High baryon fraction phase (HBFP) 

From Figures 4 and 6, we see that almost all the baryonic gas with F b > 1 is hot 
(T > 10 5 K) and moderately dense (pdm < 5). From these results, we can define a special 
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phase of baryonic gas with the indicator Fb > 1: it can be called high baryon fraction phase 
(HBFP). Baryonic gas in the HBFP is characterized by two properties: (1) it is located 
mostly in the moderately dense regions with Pdm < 5, and (2) its temperature is larger 
than 10 5 K (see §4.2). Conversely, Figure 7 indicates that almost all the gas with Pdm < 5 
and T > 10 5 K has mean Fb larger than 1. Figure 8 is another version of Figure 7, which 
shows the p£>M-dependence of the mean baryon fraction for baryonic gas with T < 10 5 K and 
T > 10 5 K. It also indicates that almost all the gas with T > 10 5 K and p DM < 5 are F b > 1. 
More precisely, 89% baryons in the T > 10 5 K and p DM < 5 regions are F b > 1. Thus, the 
phase of Fb > 1 is approximately equivalent to the thermodynamic condition T > 10 5 K and 
Pdm < 5. The HBFP is a special phase of baryonic gas formed as a result of the nonlinear 
evolution of baryon-dark matter systems. 

Figure 9 is the same as the z = panel of Figure 8, but with data smoothed on scales 
0.26, 1.04, and 4.17 h~ l Mpc. All the basic features of Figure 9 are the same as those of 
Figure 8. That is, Fb < 1 is generally in regions with high density and high temperature, 
while F b > 1 is mostly in T > 10 5 K, but with lower Pdm- In the range 0.03 < pdm < 1, the 
curves of Fb in Figure 9 are clearly less sensitive to smoothing scales. This shows that shot 
noise does not affect our conclusions about the mean baryon fractions (see Appendix B). 

On large enough scales, the cosmic clustering can still remain in the linear regime, 
and therefore one may expect that the inhomogeneity of baryon fraction distribution would 
disappear on large scales. The mean baryon fraction within a large radial region around 
a cluster has to be asymptotically approaching 1. Figure 9 shows that Fb is still quite 
inhomogeneous even at a scale as large as 4.17 h^ 1 Mpc. Therefore, the asymptotic radius 
for F b — > 1 should be larger than 4 h" 1 Mpc. 

In terms of observations, the HBFP is very different from other phases of the baryonic 
gas. As for the mass density, the HBFP is about the same as the IGM for Lja forests, 
which is given by the absorption of HI in the regions 0.5 < pou < 5 (Bi & Davidsen 1997). 
However, the HBFP cannot be detected by the Lya forests of QSO's absorption spectrum, 
because the fraction of HI is too low to be seen when temperature T is higher than 10 5 
K. As for the temperature, the HBFP is about the same as the so-called WHIM (warm hot 
IGM) (Cen & Ostriker 1999), which is generally defined as baryonic gas with temperature 
10 5 — 10 7 K, located in the regions close to the filaments, with mass density 5 < pdm < 200 
(Dave et al. 2001). The WHIM would be a source of soft X-ray background. However, the 
number density of the HBFP gas is too low for it to be a source of soft X-ray emissions, and 
hence it cannot be detected via the soft X-ray observations either. For these reasons, the 
HBFP baryons are out of the baryon budget counted with the current observations. The 
contribution of the HBFP to the total cosmic baryon budget is given in Figure 10. The mass 
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fraction of HBFP in the baryon budget is about 2.5% at z = 3 and increases to 14.4% at 
the present day (z — 0). Therefore, the HBFP would not be a negligible component of the 
missing baryons. 

5. Conclusions and Discussions 

In the nonlinear regime, a system consisting of collisional baryons and collisionless dark 
matter is generally characterized by strong shocks and discontinuities in the baryonic fluid. 
The shocks outgoing from high-density regions generally slow down the infall motion of 
baryons from low-density to high-density regions. Consequently, the baryon fraction in 
lower overdense areas is higher than the cosmic value, while in higher overdense areas it is 
lower than the cosmic value. We use iV-body/hydrodynamic simulation samples produced 
by the WENO code to quantitate the deviation of the baryon fraction from cosmic value. 
We conclude that the overall baryon matter of massive halos like clusters is lacking by about 
10%-20%. We also find that the HBFP is composed of the gas with T > 10 5 K and is located 
in moderate-density (pdm < 5) regions. 

Thus, in terms of the baryon fraction, gas in massive halos, like clusters, is in the 
phase of low baryon fraction, while gas with p DM < 5 and T > 10 5 K is in HBFP. Both 
the low baryon fraction phase and HBFP are formed by the same mechanism of the shock- 
caused separation between the baryonic gas and dark matter. About 14% baryons in the 
universe today are presumably hidden in the HBFP. The HBFP can be traced neither by 
QSO absorption spectrum nor by X-ray emissions. However, the ionized electrons in the 
HBFP would be capable of scattering CMB photons and might generate secondary cosmic 
temperature fluctuations. The Sunyaev-Zel'dovich effects might be promising for detecting 
the existence of the HBFP. 

We should point out that star formation and their feedback on the IGM evolution are 
not considered in our simulation. Generally speaking, there are two types of feedbacks: (1) 
photoionization heating by the UV emission of stars and AGNs and (2) injection of hot gas 
and energy by supernova explosions or other sources of cosmic rays. The photoionization 
heating actually can be properly considered, if the UV background is adjusted by fitting the 
simulation with the observed mean flux decrement of QSO Lya absorption spectra (Feng et 
al. 2003). The effect of injecting hot gas and energy by supernovae is localized in massive 
halos, and thus it may change some results with clusters but does not affect the IGM in low- 
and moderate-density areas. Therefore, the features of the HBFP would not be significantly 
affected even if considering the effect of star formation. 
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A. Smoothing with scaling functions 

Consider a one-dimensional density fluctuation 5(x) on a spatial range from x = to L. 
We divide the space into 2 J segments labeled 1 — 0,1, ...2 J — 1, each of size L/2 J . The index 
j is a positive integer and gives length scale L/2- 7 . The larger the value of j, the smaller 
the length scale. The index / represents the position and corresponds to the spatial range 
lL/2 j < x < (I + l)L/2 j . Hence, the space L is decomposed into cells 

The discrete wavelet is constructed such that each cell (j, I) supports a compact function, 
the scaling function <fij^(x). In our calculations, the Daubechies 4 (D4) wavelet (Daubechies 
1992) is used. The scaling function satisfies the orthogonal relation 

j <t>i,i{ x )<t>jA x ) dx = <^z'> (Al) 

where 5 K is Kronecker delta function. The scaling function 4>j,i{x) is a window function on 
scale j centered at the segment /. The normalization of the scaling function is f (fijj(x)dx = 

For a field p{x), its mean in cell (j, I) can be estimated by 

j^^(x)^x)d g = /gy/ 2 

where is called the scaling function coefficient (SFC), given by 



I p(x)4>jj(x)dx. (A3) 
A one-dimensional field p(x) can be decomposed into 

2^-1 

^) = E e ^v(^) + 0(>j)- (A4) 



1=0 
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The term 0(> j) in equation (A4) contains only the fluctuations of the field F(x) on scales 
equal to and less than L/2 J . This term does not have any contribution to the window 
sampling on scale j. Thus, for a given j, the one-point variables pjj or e P t (/ = 0, 1...2- 7 — 
1) give a complete description of the field F(x) smoothed on scale L/2K As one-point 
variables the e P t are similar to the measure given by count-in-cell technique. However, 
the orthonormality equation (Al) ensures that the set of pjj or t p - { does not cause false 
correlations. When the "fair sample hypothesis" (Peebles 1980) holds, the average over the 
ensemble of the random field can be estimated by averaging over modes (j, I). 



B. Errors of Poisson sampling 



Consider a random field p M (x) = p[l + where 5(x) = [p{x) — p]/p, and p is the 

mean. Obviously, (S(x)) = 0. The observed or simulated field p(x) is considered to be a 
Poisson sampling of the field p M (x). The characteristic function of the p(x) is 



Z [ e * /,(*)«(*)«**] = exp |y dxp M (x)[e lu{x) - 1] j , 



and the statistic of p(x) is given by 



(p( Xl )...p(x n ))p = — 



5 n Z 



5u(xi)...5u(x n ) 



u=0 



where {-..)p is the average for the Poisson sampling. We then have 

(p(x)) P = p M (x), 

and 

(p(x)p(x')) P = p M (x)p M (x') + 5 D (x - x')p M (x). 
Subjecting equations (B3) and (B4) to the transform equation (A2), we have 

(Pj,i)p = Pi 

and 

2 J ' f 

(Pj,iPj,i)p = pfipfi + J p M (x)<p j , l (x)<p jyl (x)dx. 
Therefore, the Poisson sampling error of the measurement p^i can be estimated by 

2? 



= (Pj,iPj,i)p ~ (Pj,i)p 



p (x)(j) jt i(x)(j) j) i(x)dx. 



(Bl) 



(B2) 

(B3) 
(B4) 

(B5) 
(B6) 

(B7) 
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If we use the normalized density variable, equation (B7) gives 

„, = > ez\ - /at V = %\ + * f ^ Mx)Mx)ix . ( B8) 



Since the Poisson samplings for different modes (j, I) are uncorrelated, the Poisson sampling 
error of the measurement of (1/iV) J2 N Pj,i, which is the mean of p^i over N modes, is 



N N J r 



1/2 

(B9) 



Considering that (S(x)) = 0, the second term on the right-hand side generally is negligible, 
and we then have 

^^Gz/bf . (B10) 

where N is the number of modes. The factor (L/2- J )p is the mean mass, or mean number of 
particles in the cell L/2 j . A smoothed field takes a larger (L/2 j ). Thus, if a statistical result 
is weakly dependent on the (L/2 J )-smoothing, the effect of Poisson error is negligible. 

In some calculations, we only choose the modes for which the density is restricted to 
a given range. The average over these modes may not give (S(x)) = 0. In this case, the 
S(x) term in equation (B8) is not negligible. For instance, if (5(x)) ~ 100, the error would 
be increased by a factor of 10. Nevertheless, the contribution of the 5(x) term is also 
proportional to If a statistical average over 500 modes with mean number density is 

larger than 0.03, the Poisson error will not be larger than 25%. Therefore, in Figures 7-9, 
the shot noise in the range of pdm > 0.03 are negligible. 
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Table 1: Mean Fb for Mass Density puu and Redshift z 



Pdm (pdm) z=0 z=0.5 z=1 z=2 z=3 z=4 



0.03 - 0.5 1.23 1.20 

0.5-5.0 1.17 1.11 

5.0-100 1.11 1.00 

100 - oo 0.80 0.89 



1.17 1.12 1.09 1.06 

1.07 1.05 1.04 1.03 

0.99 0.96 0.92 0.89 

0.81 0.83 0.72 N/A' 



Data are not available for p^M > 100 at z = 4. 



-17- 




Fig. I. — Example of the one-dimensional distribution of F b (x) at z=0 and 0.5, where x is 
the comoving coordinate. The corresponding temperature and density fields are also shown 
in the bottom panels. 
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Fig. 2. — Density contour plots of dark matter (right panel) and baryonic gas (left panel) 
for a slice of 0.26 hr l Mpc thickness at z — 0. The solid contours encompass the overdense 
regions with p = e i//2 , i = 0,1,2..., (p is normalized to 1). The dotted lines represent the 
boundaries of the underdense regions with p = e~^ 2 , i = 1,2.... 



as 



• 




'V*' 1 1 

1 5 v. 



Fig. 3. — Contour plots of baryon-to-dark ratio Fj, (right panel) and temperature T (left 
panel) for the same slice as in Fig. 2. The solid contours represent, respectively, the regions 
with F b = e i/3 , % = 0, 1, 2, ... and T = e l/2 x 10 5 K, % = 0, 1, 2, .... The dotted lines represent 
F fe = e"^ 3 , i = 1,2, ... and T = e"^ 2 x 10 5 K, % = 1, 2, ... regions. 
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Fig. 4. — Relation between F b (x) and Pdm(x) for each point x at redshifts z—6, 4, 2, 1, 0.5, 
and 0. In each panel, the data consist of 19,200 points randomly drawn from the simulation 
samples. The density pdm is in units of Pdm- 



-20- 



uf 4 



Lif 



L=4.17h' 1 Mpc 


L=1.04h 1 Mpc - 

H-i-t| — ' i iiw| — i i II 1 1 1 1 1 i i iiiiii| — i i Miiij — i i inn 


^ — i — i — i— »h — i — i — i — h-i— 

L=0.2Gh' 1 Mpc '. 

■ ■ ■■ ... : 


L=0.13h 1 Mpc '. 



Q I i_i i_id ■ ' ' • • I Lu_d . t . re..7TT ,V1 ^ * * 

10 3 10 2 10' 1 1 10 10 2 10 3 10 2 10' 1 1 10 10 2 10 3 
P*1 0dm 



Fig. 5. — Relation between F b and pdm for data smoothed on scales 0.13, 0.26, 1.04, and 
4.17 h- 1 Mpc at z = 0. 
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Fig. 6. — Relation between F^x) and temperature T(x) for each point x at redshifts z—A, 
2, 1, and 0. The data consist of 19,200 points randomly drawn from the simulation samples. 
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Fig. 7. — Averaged vs. temperature T. The three curves are for (1) the total samples 
(solid lines), (2) the data with p DM < 5 (dotted lines), and (3) the data with p DM > 5 
(dashed lines). 
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Fig. 8. — Averaged F b vs. Pdm- The three curves are for (1) the total samples (solid lines), 
(2) the data with temperature T < 10 5 K (dashed lines), and (3) the data with T > 10 5 K 
(dotted lines). 
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Fig. 9. — Same as Fig. 8, but only for z — 0, and the data are smoothed on scales of 0.26, 
1.04, and 4.17 hr l Mpc, respectively. 
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Fig. 10. — The z-evolution for the mass fraction of the HBFP (high baryon fraction phase) 
in the total baryonic matter. 



